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Isospin-dependent clusterization of Neutron-Star Matter 
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Because of the presence of a liquid-gas phase transition in nuclear matter, compact- 
star matter can present a region of instability against the formation of clusters. We 
investigate this phase separation in a matter composed of neutrons, protons and 
electrons, within a Skyrme-Lyon mean-field approach. Matter instability and phase 
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Q\ ' properties are characterized through the study of the free-energy curvature. The 



effect of /3-equilibrium is also analyzed in detail, and we show that the opacity to 
neutrinos has an influence on the presence of clusterized matter in finite-temperature 

o 



proto-neutron stars. 
PACS numbers: 



I. INTRODUCTION 



Nuclear matter is known to present a phase transition o: 



at sub-saturation density up to a critical temperature 



the liquid-gas type, which occurs 



aas. 



In laboratory experiments, 
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this phase diagram can be explored only through heavy-ion collisions; in particular the multi- 
fragmentation phenomenon, observed for collisions around the Fermi energy and above, has 
been interpreted as the occurrence of the nuclear liquid-gas phase transition rounded by 
finite-size effects . From the theoretical point of view, the isospin properties 

of the phase diagram have been studied extensively Iiol. H| and we know that the phase 
transition concerns isospin-symmetric (p n = p p ) as well as asymmetric (p n 7^ p p ) matter. 
Then the transition phenomenology may also affect the physics of compact stars, which are 
made of dense, neutron-rich matter [3, [3] . In particular neutron-star properties may be 



influenced by the presence of exotic "pasta phases" in their inner crust [14] . 
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Density fluctuations due to the liquid-gas instability may also play a role in the dynamics 
of supernova explosions, by affecting neutrino propagation Q] ■ 

In this paper, we study the instability of star matter against the liquid-gas phase transi- 
tion. At the densities of interest for our study, it can be considered as a medium composed 
of neutrons, protons and electrons, under the constraint of strict electroneutrality. We will 
describe nuclear matter (neutrons and protons, denoted "iVM") in a mean-field approach, 
using different Skyrme effective interactions. The Skyrme-Lyon parameterization Sly230a, 
already employed for astrophysical applications in ref. 20 L is taken as a reference, since it has 



been especially fitted to describe neutron-rich matter [2JJ, |22j • It is compared to two other 
Skyrme forces : the earlier SGII23J, whose parameters are optimized to describe nuclear 
collective motion, and SIII, one of the first and more simple Skyrme parameterizations 24]. 
Fig. [1] shows the density behavior of these different forces in the isoscalar (upper part) as 
well as isovector (lower part) channel. 

We can see that realistic interactions, fitted on measured nuclear properties as Sly230a 
and SGII, are perfectly superposed for symmetric nuclear matter, while SIII presents an ab- 
normally large incompressibility, given by the curvature at saturation : = 9pg g k 

p P=P0 

The density dependence of the symmetry energy is less well constrained by experimental 
data[10(], and shows very different behaviors with the different interactions. 

Electrons are described as a Fermi gas, whose density is set equal to proton density. 
Nuclear matter with electrons (representing star matter) will be denoted " NMe" . 

In section [TTJ we present nuclear-matter thermodynamics and compare the spinodal- 
instability region for NM and NMe. In section IHIl we consider the instability of homoge- 
neous NMe against finite-size density fluctuations, leading to matter clusterization. Some 
qualitative features of the clusters are deduced from the free-energy curvature properties of 
homogeneous matter in the framework of a spinodal-decomposition scenario, where cluster 
characteristics are determined by the dominant fluctuation mode at the onset of phase sepa- 
ration. Finally, in section [IV] we investigate the effect of neutrinos on star-matter instability 
when /3-equilibrium is assumed. 
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FIG. 1: Energy proporties at zero temperature for the different Skyrme forces used in this paper. 
Saturation points are indicated by dots. Top : energy per nucleon of symmetric matter. Bottom : 



density dependence of the asymmetry coefficient a s 
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with I = ps/p. 



II. GLOBAL INSTABILITY OF HOMOGENEOUS MATTER 

In this section, we establish the thermodynamic instability properties of infinite and 
homogeneous matter. For a fixed temperature, the free-energy curvature in the density plane 
{pni Pp) determines the spinodal region. Inside the spinodal, a thermodynamic fluctuation of 
the homogeneous system can lead to a reduction of the free energy : this fluctuation is then 
amplified, leading to phase separation. Since the resulting phases are both homogeneous and 
infinite, any even infinitesimal charge density would lead to a divergent Coulomb energy: 



25) 



the thermodynamic limit thus imposes that each phase is exactly neutral 

In the actual cooling of a proto-neutron star, the characteristic time scales may be long 
enough respect to the typical relaxation times associated with the strong interaction, for the 
system to be in thermodynamic equilibrium at all steps of its evolution. If this is the case, 
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the spinodal zone is not directly pertinent for star-matter thermodynamics. However the 
spinodal region contains points of locally inversed curvature, and is always by construction 
included inside the region of globally inversed curvature, namely the coexistence zone 
Therefore our results on the thermodynamic instabilities will keep being relevant for the 
equilibrium mixed phase. 

We will first evaluate the spinodal region in nuclear matter, in the absence of Coulomb 
interaction. Then, electroneutrality will be insured by the presence of electrons, which 
results in a strong quenching of instability. 



A. Nuclear-matter spinodal region 

Nuclear matter involves an infinite quantity of nucleons interacting through the nuclear 
force. The Coulomb interaction is not taken into account in this model, but neutrons and 
protons are distinguished by their isospin. 

In a mean-field approach, each nucleon is an independent particle evolving in a potential 
resulting from the averaged effect of other nucleons. Let us denote H the hamiltonian of the 
system, and p q the density matrix of each nucleon species, the index q representing neutrons 
(q = n) or protons (q — p). The mean field W q is defined by the relation 5(H) = Tr(W q 5p q ). 
Within this approximation, the average energy density of homogeneous nuclear matter H h = 
(H) /V is a functional of one-body densities only: particle densities p q and kinetic densities 
r q = (p 2 ) q /h 2 . It is useful to introduce isoscalar and isovector densities : 

P = Pn+Pp, T = T n + T p , 
P3 = Pn~ Pp , T 3 = T n - T p . 

Using a Skyrme effective interaction, the energy density of homogeneous, spin-saturated 
matter with no Coulomb interaction, can be written as : 

H = K + 7io + H3 + H e ff (2) 
In this expression, fC is the kinetic-energy term, TLq a density-independent two-body term, 
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Ti 3 a density- dependent term, and Heff a momentum- dependent term: 

/V - 2m ' 

Ho = C p 2 + A)P§ 

n 3 = c 3 p°+ 2 + D 3P °pi 

7-teff = Ceffpr + D effP3 T 3 
where the coefficients Cj an Di are linear combinations of the standard Skyrme parameters 
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211 ] . The mean field is given by the expression : 



where m* is the effective mass of particles of type q. The individual energy levels read 

„2 

e q = 2m* Uq- m the grand-canonical formalism, they are occupied according to Fermi- 
Dirac statistics with occupation numbers: 

n q (p) = 1/ [1 + exp((3(p 2 /2m* q + U q - //,))] , (5) 

where /3 is linked to the temperature by (3 = 1/T, and p, q is the chemical potential of 
particles of type q. All one-body densities are then given by Fermi integrals and the mean- 
field grand-canonical partition sum can be expressed in the same way as for an ideal gas: 

— = 2 2 r J Q Ml + e )— dp = ^—{fr qt (6) 

where //' = \i q — U q . 

For fixed values of p n and p p at a given temperature, the value of chemical potentials 
as well as kinetic densities are obtained by iteratively solving the self-consistent relations 
between p q , pL q and m*. We can then calculate the entropy as a function of densities. Indeed 
the mean-field approach approximates the real entropy S by the corresponding mean-field 
entropy S° 2(| : 

8 = S/V ~ S /V = \nZ /V + P({W) Q /V - VnPn - P P P P ). (7) 

For a given temperature T = statistical equilibrium corresponds to the minimization 
of the free energy F = Vf, related to the entropy S by a Legendre transform. Per unit 
volume, this relation is / = 7i h — s//3 and / is a function of nucleon densities p n and p v . If the 
free-energy surface presents a concavity, it can be minimized by linear interpolation (phase 
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mixing). The equilibrium configuration corresponds to a phase coexistence determined by 
Gibbs contruction, which is a several-dimension generalization of the well-known Maxwell 
construction {<]]. 

We now turn to spinodal instabiliy. Let us consider an homogeneous system submitted 
to a thermodynamic fluctuation, consisting in a separation into two portions distinguished 
by an infinitesimal density shift. Both portions correspond to infinite and homogeneous 
systems. If the fluctuation leads to a reduction of the global free energy, it will be amplified 
in the successive time evolution: the initial system is locally unstable. The spinodal region is 
then defined as the ensemble of unstable points, i. e. points for which the free-energy surface 
is convex in at least one direction. This is determined by the curvature matrix : 



This matrix is defined for each point (p n ,p p ), and its lower eigen- value C< gives the mini- 
mal free-energy curvature in this point. The eigen- vector u< associated with the minimal 
curvature, named "instability direction", gives the direction of phase separation. 

We have represented on the left part of figf2]the minimal curvature C < as a function of 
total density, for different values of the proton fraction Z/A, at a fixed temperature of 10 
MeV : some of the curves have a negative region. On the right part, the minimal value 
of C < obtained for each asymmetry is represented as a function of Z/A. This shows that, 
for a finite range of asymmetry, homogeneous matter presents a spinodal instability, i.e. a 
negative curvature of the free energy. This is observed until a critical temperature of about 



Let us now take the electric charge into account. Matter must be strictly neutral at a 
macroscopic scale, in order to avoid the divergence of Coulomb energy. In this paper, we 
consider star matter where the proton charge is neutralized by a Fermi gas of electrons. 
Studying thermodynamic fluctuations, each portion in which the system separates is sub- 
mitted to the constraint of strict electroneutrality : the electron density is then everywhere 
exactly equal to the proton density, p e = p p . In this case, electrons do not add any degree 





14 MeV. 



B. Electron quenching of instability 
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FIG. 2: Minimal free-energy curvature obtained with the Sly230a interaction. Ordinary nuclear 
matter (iVM)presents a spinodal region, which is suppressed when an electron gas is added to 
neutralize the proton charge (NMe). Left : lower eigen- value as a function of the total density p, 
for different Z/A from Z/A = 0.5 to Z/A = 0, equally spaced of 0.1. Right : minimal curvature as 
a function of Z/A. Dotted line : MNe at T = 0. 



of freedom, and the free-energy curvature matrix remains two-dimensional : 

nh I dfl n /dp n dp n /dp p \ 

^ 'NMe ~ I I V*) 

y d/ip/ dp n d/j p / dp p + dp e / dp e J 

The presence of electrons contributes to this matrix via the electron susceptibility Xe 1 — 
dp e /dp e . This term is readily evaluated in the approximation of a degenerated, ultrarela- 
tivistic Fermi gas. Because of the low electron mass, these approximations are valid until 
very low densities. For a degenerated Fermi gas, the Fermi momentum pp e is given by the 
relation : 

Pe = Pp= 3^W" (10) 
In the ultra-relativistic approximation, the chemical potential is equal to the Fermi momen- 
tum p e = Pf b c, so it is related to the density by : 

p e = /ic(37r 2 p e ) 1/3 , (11) 

which gives the electron contribution to the curvature matrix : 

dp e /d Pe = Xe' = H^- 2 ) 1/3 - (12) 
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In the density region of interest for our study, the numerical value of Xe 1 * s large i n 
comparison with the other terms of the free-energy curvature matrix 35(. As a result, 
adding an electron gas leads to a drastic effect on the eigen-modes of the matrix. Because 
of its positive sign, it leads to a quenching of the instability. This is illustrated on figfS] : 
the middle part presents the minimal curvature as a function of total density for MNe with 
different values of Z/A, at T = lOMeV. The lowest minimal curvatures are reported on 
the right part as a function of Z/A. No negative value of the curvature is found, and the 
same is true at zero temperature. Since the extension of the spinodal region monotonically 
increases with decreasing temperatures, this result implies that NMe does not present any 
thermodynamic instability. 

The value of (7™ m depends on the details of the nuclear interaction used for the calculation 
of neutron and proton chemical potentials. For one given point, the sign of C < is given by 
det C^Me = det C/vm + ® p-nl^n * d Pe fJ, e . It can be negative if det C^ M has a sufficiently strong 
negative value, ie nuclear-matter thermodynamic instability is strong enough to balance 
the positive term due to the electron gas. This is does not happen for any of the three 
Skyrme forces we have considered, as shown by figure [3]: the instability is totally suppressed 
by electron incompressibility. However, this result is model dependent. For example, in 
reference 27| the spinodal instability is studied in a relativistic-mean-field approach, using 
the NL3 set of parameters : the authors obtain a small thermodynamic instability region 
in the case of nuclear matter with electrons. This region, located at low density, is very 
reduced with respect to nuclear-matter spinodal, as expected from electron quenching. 

III. CLUSTERIZATION OF NMe 

In the previous section we have shown that the electron gas in stellar matter, because of its 
high incompressibility, acts against the formation of macroscopic charge dishomogeneities, 
thus completely quenching the liquid-gas phase transition, at least in the mean-field approx- 
imation adopted here. This does not mean that we should expect stellar matter at subsatu- 



28| that 



ration density to be homogeneous. Indeed it has been recognized since a long time 
the ground state of neutron-rich stellar matter at low density should be a Coulomb lattice 
of exotic nuclei embedded in a neutron gas, meaning that we should expect a solid-liquid 



transition between the core and the inner crust of a neutron star 1 3(31 . Numerical simulations 
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FIG. 3: Suppression of spinodal instability in nuclear matter with electrons, for different Skyrme 
forces at T = 0. Top : minimal free-energy curvature obtained for a given proton fraction, as a 
function of Z/A. Bottom : density at which this minimum is obtained. 



performed with (semi)-classical models [ly, ll7| indicate that such dishomogeneous phases 
should persist at finite temperature and in a large range of asymmetry. Clusterisation of hot 
stellar matter is also considered in a statistical approach in ref. 29j, |30[]. A way to address 



the possible existence of a phase of clusterized matter in mean-field approaches, is to study 



the instability of nuclear matter with electrons against finite-size density fluctuations 



31|. If 



such fluctuations lead to a reduction of the free-energy, they will be amplified through the 



process of spinodal decomposition 



321 ] until the formation of structures of alternatively high 



and low density. Homogeneous matter is then unstable against clusterization. Clusters cor- 
respond to the high-density regions, which can be seen as liquid nuclear drops immersed in a 
nuclear gas. Considering the short time scales associated with spinodal decomposition, the 
cluster characteristics are determined by the free-energy curvature properties at the onset 



of instability, when density fluctuations have small amplitudes 



32|]. 
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A. Free-energy variation with finite-size fluctuations 

The instability region of homogeneous NMe against clusterization is determined by 
studying how the system free energy changes, when a finite-size density fluctuation is intro- 
duced. 

If fluctuations occur on a finite microscopic scale, electron and proton density can fluctu- 
ate independently : only their mean values are constrained to be equal, insuring macroscopic 
electroneutrality. Fluctuations thus affect independently the three species of the medium 
(neutrons, protons and electrons), whose densities become : 

Pq = Pq + $Pq> (!3) 

with q G {n,p, e}. Each density variation can be expressed by a Fourier transform 

5p q = [ dka q {k)j tp , (14) 



— * — * — * 

with a q (k) = a*(—k) to insure that Sp q is real. Since the different wave vectors k are 
decoupled in the global free-energy variation, the problem reduces to the study of plane- 
wave density fluctuations: 

5p q = A q e ik ' r + c.c, (15) 

where each species is associated with a different amplitude. Thermodynamic fluctuations, 
studied in the previous section, can be obtained as the limit k — > of this model. 

To evaluate the free-energy variation, we consider a Thomas- Fermi approximation, i.e. 
the density variation is supposed smooth enough to allow in each point the definition of 
a Fermi sphere corresponding to the local density. Then, in each point of density p q (r) = 
Pq + Spq-, the local bulk term of the free energy is equal to the free energy f h of an infinite 
homogeneous system at the same density. The global bulk free energy of the system is the 
space average of this local term : 

f=\\ f h ({p q (r)})dr. (16) 
In the small-amplitude limit, the integration leads to f b = f h ({p q }) + 5f b , with : 
, F i._S- A ' A 'i + A - A ' ( 92 f" \ ^ A.A' + A-A, /d n \ 
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where i,j G {n,p, e}. First-order terms have vanished in the integration because the average 
density variation is zero. 

The variation of the entropy is contained in the bulk term, since entropy depends only 
on the local density |26|]. However, in the case of a finite wave number k, the energy density 



is modified by two additional terms 3jj, arising from the density-gradient dependence of 
the nuclear force, and from the Coulomb interaction. Denoting these two contributions <5£ v 
and 5S C respectively, the free-energy variation is : 

Sf = Sf b + S£ v + S£ C . (18) 

Let us now explicit the gradient and Coulomb contributions. In the presence of density 
gradients the nuclear energy density has the form 7i = T~C h + ?^ v , where 7i h is given by 
eq.(j2]) and the density-gradient term TC 7 is expressed as: 

H v = Cj n {V Pn f + Cj p (V Pp ) 2 + 2Cl V Pn ■ V Pp . (19) 



Coefficients are combinations of the Skyrme parameters given for example in ref. 2lfl : 



Cl = Cl = £ [ tl (i _ Xl ) - t 2 (l + x 2 )) , 
CJ = C p v n = i [3^(2 + x x ) -t 2 (2 + x 2 )} . 

The global contribution of this term to the energy density is given by the space average: 

8£ v = i J H v (7)dr= k 2 (AiA* + A* A,) Cj, (21) 

i,3 

where i, j G {n,p}. 

The Coulomb contribution 5H C is due to the independent density fluctuations of charged 
particles {e,p}. Denoting % = y/Aneoei the electric charge of a particle of type i, this creates 
the charge distribution : 

Pc(?) = J2 e * A ^+ c - c - ( 22 ) 



The consequent Coulomb energy per unit volume is : 

i f pJiikM dtldt2 = t 

2V J |r 2 — ri| k 2 

where i, j G {e,p}. 
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Summing the contributions ([IT]) . ([21]) and (|23l) . we have the total free-energy variation, 
expressed to the second order in A q : 



A ■ A* -L A* A ■ 



d Pj 



47re,e 



{p° q } 



k 2 



(24) 



which can be written in a matrix form in the three-dimensional space of density fluctuations 

^4 = (A n , A p , A e ) : 

5f = A*C f A, 



(25) 



where 



C* 
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(26) 



is the free-energy curvature matrix. The bulk term of the curvature is conserved in this 
new expression, with the difference that the electron term now belongs to an independent 
space. The density-gradient part of the nuclear interaction adds a term proportional to k 2 , 
while the Coulomb interaction induced by the plane-wave charge distribution adds a term 
inversely proportional to k 2 . 



B. fc-dependent instability 

An instability against matter clusterization corresponds to a negative free-energy curva- 
ture in the space of density fluctuations. It is studied through the /^-dependent curvature 
matrix . The eigen-mode of minimal curvature is represented in figJUfor one given point of 
average density (p°,Pp), as a function of k, by three quantities : its eigen- value, its direction 
in the (p p , p n )-plane and in the (p p , p e ) -plane. 

The upper part of figJH gives the minimal eigen- value C< of the curvature matrix. We 
observe that, for the chosen average- density point, there is an interval of k allowing a negative 
free-energy curvature in the space of density fluctuations : the corresponding plane wave 
will then be amplified leading to matter clusterization. The lower limit of the instability 
interval is due to the Coulomb term (oc 1/A: 2 ), which tends to suppress long wavelenght 
instabilities which would lead to the formation of a macroscopic charge. The upper limit 
is determined by the nuclear density-gradient term (oc k 2 ), which acts against the surface 
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FIG. 4: Eigen-mode of the minimal free-energy curvature for p = 0.05/m~ 3 , Z/A = 0.3, T = 5 
MeV, as a function of the wave number k, calculated with the Sly230a Skyrme interaction. Top 
: eigen-value. Middle : associated eigen-vector in the nuclear-density plane. The iVMe-curves 
are compared to two limiting cases corresponding to a zero and infinite incompressibility of the 
negatively charged gas. When the curvature is negative (full lines) this eigen-vector gives the 
phase-separation direction. The dots give the points of zero curvature. The dashed line gives p%/ ' p 
(direction of constant Z/A). Bottom : same as the middle part, in the (p p , p e )-plane. The dashed 
line gives do^/Sor. when the nuclear force is zero. 
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formation associated with short wavelenght fluctuations. Thus, both effects arising from 
finite-size fluctuations lead to a reduction of the instability, with respect to the case of 
thermodynamic phase separation in uncharged nuclear-matter. 

Let us now consider the phase-separation direction, which is given by the eigen-vector 
associated with the minimal curvature. The middle part of figJH represents this direction 
projected in the (p n , p p )-plane. Comparing this direction to the line of constant proton 
fraction, we can see that the amplification of fluctuations will lead to a dense phase (clusters) 
which is more symmetric than the surrounding gas. This is the well-known phenomenon 
of isospin distillation 9], resulting from the minimization of the nuclear symmetry energy. 
The lower part of the figure represents the eigen-vector direction in the (p p , p e ) plane, i.e. it 
represents the mobility of electrons. This behavior is essentially determined by the relative 
weight of coulomb interaction and electron incompressibility, as shown by the fact that the 
complete neglect of nuclear interaction (dashed line) changes the results only slightly. At 
low A;-values, this quantity is dominated by Coulomb energy, which forces the electrons to 
closely follow the protons in order to lower the electric charge. Conversely, for high values 
of k, electron incompressibility dominates over Coulomb interaction, and the electron gas 
tends to stay a uniform background. 

In order to better quantify the effect of electrons, the phase-separation direction is com- 
pared to two extreme idealized cases represented on the same figure. In both cases, the 
proton charge is neutralized by a gas of negative charges : the incompressibility of this gas 
is set to either zero or infinity. In the first case, the Coulomb energy is exactly cancelled : 
the resulting curves show the effect of nuclear density-gradient term, the point at k — cor- 
responding to homogeneous nuclear matter as discussed in section [TTAl In the second case, 
the negative gas forms a strictly homogeneous background, which suppresses the monopolar 
term of Coulomb interaction. In reality the electron gas has a finite incompressibility, related 
to the susceptibility x" 1 . As noticed in section HTT?| the high value of x~ x strongly suppresses 
the electron degree of freedom. As a result, the electron dimension in the curvature matrix 
( |26l) has a only a perturbative effect with respect to the case of a uniform negative back- 
ground. Since the (small) electron mobility allows to lower the Coulomb energy, it makes 
the fluctuation more favorable and leads to a slight increase of the instability interval. 

Fixing different values of k, we have represented in figj5] the corresponding instability 
regions (fc-spinodals) in the density plane (p n , p p ) at a given temperature. Each point of 
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this plane gives the initial densities of the homogeneous system, from which fluctuations 
start. The effect of the finite electron incompressibility is more visible at lower k, where the 
Coulomb term dominates. For a fixed k, the perturbation increases as \~ l decreases, i.e. 
towards higher values of p° e = p®. Starting from k ~ 80 MeV/c the electron gas can safely 
be represented as a fixed background. The curves presented on figfS] can be qualitatively 



compared with the results of ref.[27| (obtained in a RMF approach at zero temperature). 
Nuclear-matter spinodal region are very similar in both cases. When finite-size instabilities 
are considered, a similar behaviour is obtain for the higher values of k (k > 100). For the 
lower values, Coulomb contribution becomes dominant and is found to affect differently the 
shape of the /c-spinodals : figj5] shows a stronger reduction of instability in the proton-rich 
part of the plan, while the curves remain essentially symmetric in the RMF calculation. 
This signs a different evolution of the derivative d Pn fi n through the density plane (indeed, 
when Coulomb interaction dominates, the direction of minimal free-energy curvature gets 
closer to the neutron-axis direction). 

The /c-interval, for which NMe presents an instability in at least one density point, is 
represented on figE] as a function of the temperature T. It is maximal at zero temperature, 
and goes to zero at a limiting temperature which is about YlMeV for the realistic forces 
(almost lQMeV for SIII). 

The global instability region (at a given temperature) is defined as the envelope of the dif- 
ferent fc-spinodals : it contains all the points presenting an instability against clusterization 
for at least one value of k. It is represented in figd for three different Skyrme parame- 
terizations : Sly230a, SGII and SIII. We can see that the two realistic forces give similar 
envelopes. SIII presents an atypical behaviour, with a much wider instability region. This is 
essentilly due to the too high value of (nuclear incompressibility at saturation density) 
(see Fig.l). 

C. Influence on cluster properties 

In a spinodal-decomposition scenario, phase properties are mostly determined by the 



conditions at the onset of instability 32]. In particular matter clusterization is dominated 
by the development of the most unstable mode, i.e. the wave number giving rise to the fastest 
amplification of the fluctuation. The temporal evolution of a fluctuation of wave-number k 
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FIG. 5: fc-spinodals for different values of k at T = 5MeV. The external dotted curve is the 
thermodynamic spinodal region of nuclear matter, it is the same in all panels. The fc-spinodals 
of NMe are compared to two limiting cases (zero and infinite incompressibility of the negatively 
charged gas). The thermodynamic spinodal (dotted curve) corresponds to the case x _1 = (dashed 
curve) at k = 0. 



17 




T (MeV) 



FIG. 6: Wave-number k interval for which the density fluctuation is amplified in a region of the 
density plane, as a function of temperature. Top : NMe result (Xe 1 ^ compared to two limiting 
cases (zero and infinite incompressibility of the negatively charged gas) for Sly230a. The curve 
correponding to \~ l = goes to zero at the critical temperature of nuclear matter T u . For the 
physical NMe case, the instability interval disappears for a limit temperature Tu m < T u . This 
temperature is slightly higher for NMe. Bottom : NMe curves for three Skyrme parameterizations. 



is given by: 

5p q = A q e^-^ + c.c, (27) 

where unstable modes are characterized by an imaginary frequence. Because of the expo- 
nential amplification of such modes, clusterization is dominated by the wave number ko 
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FIG. 7: Envelope of instability against finite-size density fluctuations at T = 5MeV. For ev- 
ery point inside this region, at least one fluctuation mode can be amplified. Left : some NMe 
/c-spinodals (k values given in MeV/c) and the instability envelope (thick curve). NM thermo- 
dynamic spinodal is also given (external thin curve). Right : NMe instability envelope for three 
different Skyrme forces. 



corresponding to the highest imaginary frequence \u (ko) \ = max. The dispersion relation 
u(k) can be established in a hydrodynamic description. This is easily derived 32| in the 
case of one single fluid, from the Euler and continuity equations: 

d{mpv)/dt = - VP, (28) 
dp/dt + V(pv) = 0, (29) 

where P is the system's pressure, related to the free-energy density / by P = pdf /dp. These 
equations in the limit of small amplitude fluctuations p = p$ + 5p can be expressed as a 
diffusion equation 



d 2 p 



p CV 2 p, (30) 



dt 2 

where C = d 2 f/dp 2 \ Po is the curvature of the homogeneous- system free-energy density. A 
Fourier representation of this equation leads to the dispersion relation uj 2 = poCk 2 . When 
this curvature is negative, the frequency is imaginary : it can be expressed as u = i/r, where 
r is the time constant characterizing the development of fluctuations. 
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In our study, the instability direction for a given wave-number k in each point of the 
density plane, is associated with the eigen- value C<(fc) of the curvature matrix C. Since the 
temporal evolution of this minimal-curvature mode is given by r 2 oc l/(— C < k 2 ), the most 
unstable mode for each temperature and density can be determined by the condition —Ck\ = 
max. In the diabatic limit that the successive dynamics first approximation, be 



neglected 



321 ]. we can link the free-entropy curvature matrix of homogeneous matter, to final 



cluster properties as issued from the subsequent spinodal decomposition. In particular, an 
order of magnitude for cluster sizes is given by the half wavelength Ao/2 = 7r/k ; cluster 
isotopic composition is related to the direction of the eigen- vector w<(£;o) associated with 
C<(fco). 

The evolution of the most unstable mode through the density plane is given by figJHl for 
a fixed temperature. The upper part presents quantities linked to cluster size : the most 
unstable wave-number k , and its associated half wavelength Ao/2. For a wide range of 
asymmetries, the typical cluster size depends essentially on the total density, rising towards 
the instability boundary. Level lines circle a minimum at low density, quasi-symmetric 
matter. At these low densities the very high frequency of the most unstable mode gives 
typical sizes of a few fm, corresponding to clusters of ~ 10 particles. This result should be 
taken with care: the Thomas-Fermi approximation inherent in our treatment breaks down 
for such small clusters, and we may be facing the limit of validity of the mean field approach. 

The lower part of the figure presents properties affecting the isotopic composition of 
clusters: the direction of w<(fco) in the nuclear-density plane is given by the rapport Sp^/Sp. 
It is compared to the direction of constant proton fraction by the difference p%/p — 5p^/Sp. 
We can see that clusterization almost systematically involves portions of matter which are 
more symmetric than the global system, which is again a manifestation of isospin distillation. 

A representation of the most unstable mode is given by figJH and [10J They show some 
characteristic evolutions for the half wavelentgh and the phase separation direction as a 
function of density and proton fraction, at different temperatures and using different Skyrme 
parametrizations. We can see that we can expect clusterized stellar matter with similar 
characteristics in a very wide range of temperatures, densities and asymmetries, and these 
qualitative results appear robust respect to the chosen Skyrme parametrization. In particular 
realistic forces, irrespective of their different symmetry energy at subsaturation densities, 
lead to very close results for all observables. Increasing the temperature has the effect of 
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FIG. 8: Most unstable mode properties through the density plane, represented in level lines. 
Top-left : most unstable wave number ko(MeV/c), levels regularly spaced of 2QMeV/c. Top- 
right : corresponding half wavelength Ao/2(/m), levels regularly spaced of lfm. Bottom- left : 
most-unstable mode direction in nuclear-density plane 5p3/5p(ko), levels regularly spaced of 0.1. 
Bottom-right : difference respect to the constant Z/A direction, p^/p — 5p3/5p(ko), levels regularly 
spaced of 0.05. 
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increasing the cluster size and decreasing the distillation. The electron participation to 
phase separation is given by Sp e /Sp p : this ratio is low, and correlated with A (the larger 
clusters lead to a higher polarization of the electron gas). 

In the absence of Coulomb interaction (as for the study of normal nuclear matter), isospin 
symmetry imposes along the p n = p p axis the relation 5p 3 /<5p = p^/p = 0. In the present 
case, the presence of a proton charge hinders p p fluctuations, leading to an enhancement of 
5p 3 /Sp. Level-lines Sp 3 /5p = and Ps/p — Sp 3 /Sp = are respectively shifted to the proton- 
rich and neutron-rich part of the density plane, delimiting a region of abnormal isospin 
distillation , shown in figHH for a temperature T = 5 MeV. On the neutron-rich side of this 
region, phase separation leads to a dense phase more asymmetric than the global system, 
and a proton-rich gas. On the proton-rich side, for a large enough final density the dense 
phase is neutron rich. However, the region concerned by this anomaly is thin : for the major 
part of the instability region, the isospin distillation phenomenon of normal nuclear matter 
is not qualitatively modified by the presence of electric charge. 



We now study the case of a /3-equilibrated star matter, composed of neutrons, protons, 
electrons and neutrinos v. For simplicity, leptons will be treated in the degenerated-gas 
approximation. If this equilibrium is reached at a given temperature and density, it imposes 
a given isotopic composition, determined by the rate equality between reactions of electronic 
capture and neutron disintegration : 



When the core of a massive star collapses, the increase of density displaces /3-equilibrium 
towards a more neutron-rich composition. 

/3-equilibrium corresponds to free-energy minimization. In the presence of different species 
% of mass rrii, the total free energy density is given by : 



IV. /3-EQUILIBRIUM IN NMe 



p + e 



n + v 



(31) 



n — > p + e + v. 




(32) 



./ 



'* derivatives define the chemical potentials p\ : 



A = d Pif = d P if + rrii = pi + m, 



(33) 
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FIG. 9: Cluster properties as indicated by the most unstable mode, along two axis : one at constant 
density (left side), the other at constant Z/A (right side). Three temperatures are considered with 
Sly230a. 
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FIG. 10: Cluster properties as indicated by the most unstable mode, along two axis : one at 
constant density (left side), the other at constant Z/A (right side). Three different Skyrme forces 
are considered at T = 5MeV. 
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Pn (f m" 3 ) 

FIG. 11: Coexistence region calculated with the Sly230a interaction at a temperature T = 5 MeV. 
The grey area gives the region of abnormal distillation due to the Coulomb interaction, situated 
between the two lines 5pz/5p = and Sp^/Sp = pz/p. (dp^/Sp is the phase-separation direction 
given by the most unstable mode). 



so that /9-equilibrium corresponds to the relation 



Pp P"e 



H p + m p + /i e + m e = /i„ + m n + /!„. 



(34) 



Eq. (l34")) supposes that neutrinos participate in /^-equilibrium. This is true in the ini- 
tial post-bounce supernova dynamics, when matter is dense enough for neutrinos to be 
completely trapped in the star. In the later phase of the cooling process, matter becomes 



transparent to neutrinos 



33, 



34| , which then do not equilibrate their chemical potential with 



the one of the other species. The transition between these two regimes is smooth, and in 
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the following we will study post-bounce supernova matter considering various rates of neu- 
trino trapping. The detailed relation between density, temperature and neutrino trapping 
in supernova matter is not a settled issue, especially concerning the crust region which is 
of interest for the liquid-gas phase transition. For this reason, for our study we treat the 
trapping fraction as a free parameter, denoted x v : it is defined as the total fraction of 
neutrinos produced by electron-capture reactions which is retained in the star. We suppose 
an initial proton fraction Z/A corresponding to 56 Fe (which is the main constituent of the 
star core before collapse). For each point (p, Z/A) of the nuclear density plane, x v and Z/A 
determine the neutrino density inside the matter according to : 



Leaving free parameter, we can estimate the effect of neutrino trapping on the 

thermodynamic properties of the system. The neutrino density is linked to the chemical 
potential by the ultra-relativistic degenerate Fermi gas expression : 



Eqs. . (1351) . ( |36i) describe a trajectory in the plane of chemical potentials, which fixes the 
isotopic composition of matter at each temperature as a function of the x v parameter, linked 
to neutrino opacity. The /3-equilibrium compositions corresponding to different values of 
x v are represented in fig {12] in the space of chemical potentials (Pn, pt e — pi + pi). The 
diagonal in this space gives the /3-equilibrium curve for a matter completely transparent to 
neutrinos. Any point of the plane is /3-equilibrated for a value of p v (i.e. a number of trapped 
neutrinos) corresponding to the vertical shift respect to this line. The lines corresponding 
to x v = and 1 give the region which is in principle accessible to /3-equilibrium : this 
crosses the instability envelope, which is reported as a thick line in the chemical-potential 
plane. We can see that the matter composition at /^-equilibrium is very sensitive to the 
trapping fraction x u : only a small fraction of trapped neutrinos is sufficient to sensibly 
increase the stellar proton fraction. At zero temperature, the instability envelope is crossed 
for any value of x v . At higher temperatures of the order of T = lOMeV, relevant for the 
initial dynamical steps of the proto-neutron-star, matter at /9-equilibrium is too asymmetric 
to be unstable against clusterization if all neutrinos escape from the cooling star, while a 
fraction as small as x v = 0.01 is sufficient to lead the system to explore the instability 



Pv = x v p [{Z/A)o - Z/A] . 



(35) 




(36) 



26 




860 880 900 920 940 



860 880 900 920 940 

Mi 



860 880 900 920 940 960 



M-n ( MeV ) 



FIG. 12: /3-equilibrium in chemical-potential plane (//^,//p e = Mp + /-4)> f° r different temperatures 
and neutrino opacities x v = 0, 0.01, 0.1 and 1. The lines x u = and 1 delimit the region which is 
in principle accessible to (5 equilibrium. The NMe instability envelope is also represented. 



region. This result depends very slightly on the different Skyrme parametrizations, meaning 
that we can reasonably expect that matter will be clusterized in the whole evolution of the 
proto-neutron-star. 

A representation of /5-equilibrium composition in the plane (p,Z/A) is given by figJT3l 
for different temperatures and Skyrme interactions. Each path corresponding to a fixed 
value of x v shows a quasi-constant proton fraction. The proton fraction shows very little 
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evolution until low densities : then it first slightly decreases before going abruptly to 1 in 
the limit of zero density. Indeed, at vanishing density the mass balance m n — m p — m e > 
dominates /3-equilibrium, imposing a full decay of neutrons into protons. However, since this 
rising happens at very low density, figUS] shows that for a fixed x v /^-equilibrium essentially 
corresponds to a constant-Z/A line. 

We should stress that in this calculation nuclear matter is always considered as homo- 
geneous. When clusters are formed, /3-equilibrium does not correspond any more to the 
simple mean-field trajectories represented in figs{12] and [13j Thus the /3-equilibrium curves 
displayed by these figures represent the initial condition for matter clusterization. The insta- 
bility envelope we have represented corresponds to the region where clusters can be formed 
by spinodal decomposition. This gives a lower-limit estimate for the density interval cor- 
responding to clusterization. Indeed the region for which completely equilibrated matter 
is clusterized has to be wider than the spinodal-instability envelope. However, outside this 
envelope, matter can in principle be homogeneous and /3-equilibrated since in a situation of 
meta-stability the time needed for cluster formation can be longer than /5-equilibrium time 
scale. Let us also consider the case the case of a thermally equilibrated initially homogenous 
dense matter. If it is submitted to a fast dilatation, Z/A remains constant while density 
decreases. Then, the proportion of neutrinos trapped in the initial state, which fixes the 
constant composition Z/A, determines wether the instability envelope is crossed or not. 

Figure [14] gives a quantitative view of the neutrino fraction needed for /3-equilibrated 
matter to reach the region of instability against clusterization, as a function of temperature. 
The value x™ tn corresponds to the path at constant x v which is tangential to the instability 
envelope. The curves obtained show a slow evolution with very small values of neutrino 
fraction at low temperature, and a sharp enhancement as the limiting temperature for 
instability is approached. This profile is due to the typical evolution of the instability 
region, since, as we have seen, /3-equilibrium paths depend only weakly on the temperature. 

V. SUMMARY AND CONCLUSION 

The nuclear-matter liquid-gas phase transition results in a region of instability against 
cluster formation in low-density finite-temperature star matter composed of neutrons, pro- 
tons and electrons. This has been studied analyzing homogeneous-matter instabilities 
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FIG. 13: /3-equilibrium in the plane (p, Z/A) for different temperatures and neutrino opacities 
x v = 0, 0.01, 0.1 and 1. MNe instability envelope is also represented. 

against a plane- wave density fluctuation, in the small- amplitude limit. In a spinodal de- 
composition scenario, cluster properties as sizes and chemical compositions can be predicted 
from the characteristics of the most unstable mode at the onset of fluctuation amplification 
: they are determined studying the free-energy curvature matrix in the space of density 
fluctuations, as a function of the wave number k. 

For a cold neutron star, thermodynamic equilibrium is reached. The liquid-gas transition 
concerns the inner crust, where pasta phases are expected. Our approach does not include 
the study of such structures, since it is based on homogeneous-matter properties : however, 
the instability envelope which we obtain gives the minimal region where equilibrated matter 
is clusterized. 
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FIG. 14: Temperature dependance of the minimum neutrino opacity x™ in leading to the instability 
region on a path of /^-equilibrium. 



The out-of-equilibrium phase separation may concern proto-neutron stars issuing from 
core-collapse in supernova explosions, since the shock-wave causes abrupt changes in the 
density. At this time, star matter is hot and contains trapped neutrinos, whose propagation 



is a key issue for the understanding of supernova dynamics [18j. It has been pointed out 17] 
that neutrino in-medium propagation is affected by the eventual presence of clusters. On 
the other hand, we have found that presence of neutrinos affects the star-matter composition 
in a significant way. Indeed neutrino trapping rules the crossing of the instability region, 
which in turn leads to the formation of clusters. This interplay is sensitive to the different 
time scales implied in the explosion : temperature and density evolution, as well as the 
realization of f3 and thermodynamic equilibrium. 
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